Quantum stability of Mott-insulator states of ultracold atoms in optical resonators 
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We investigate a paradigm example of cavity quantum electrodynamics with many body systems: 
an ultracold atomic gas inside a pumped optical resonator, confined by the mechanical potential 
emerging from the cavity-field spatial mode structure. When the optical potential is sufficiently 
deep, the atomic gas is in the Mott-insulator state as in open space. Inside the cavity, however, the 
potential depends on the atomic distribution, which determines the refractive index of the medium, 
thus altering the intracavity-field amplitude. We derive the effective Bose-Hubbard model describing 
the physics of the system in one dimension and study the crossover between the superfiuid - Mott- 
insulator quantum states. We predict the existence of overlapping stability regions corresponding to 
competing insulator-like states. Bistable behavior, controlled by the pump intensity, is encountered 
in the vicinity of the shifted cavity resonance. 

PACS numbers: 03.75.Hh,05.30. Jp,32.80.Qk,42.50.Vk 



I. INTRODUCTION 

Cavity Quantum Electrodynamics (CQED) [l|, Q has 
been a key area of quantum optics since its early days 
of optical instabilities, such as optical bistability 3], till 
most recent CQED with single atoms interacting with 
single or few photons 6, 7,^,^]. 

In recent years considerable attention has been paid 
to a new regime of CQED, which we term CQED with 
many-body systems. These studies focus onto the me- 
chanical effects of the resonator field on the atomic mo- 
tion, and on the non-linearity arising from the interde- 
pendence between the cavity field and the atoms dynam- 
ics. Following the theoretical prediction of Ref. [l^], sig- 
natures of self-organization have been measured in the 
light scattered by laser-cooled atoms in a transversally- 
pumped cavity [ll|. These structures and their proper- 
ties have been theoretically studied in detail in Refs. [H, 
[l3| . In different setups, Bragg scattering of atomic struc- 
tures inside optical resonators has been experimentally 
investigated in Ref. [3]. 

While in all the cases mentioned so far the atomic mo- 
tion is essentially classical, the stability and properties of 
these structures in the quantum regime are still largely 
unexplored. This question acquires a special relevance 
in view of the recent experimental progress of CQED 
with ultracold atoms. In fact, strong atom-field coupling 
between Bose-Einstein condensed atoms (BEC) and the 
mode of a high-finesse optical cavity has been realized 
in the experiments reported in Refs. [H, [H, [O, fisj . 
Moreover, CQED techniques were used to measure pair 
correlations in the atom laser [l9j . and have been pro- 
posed for characterizing quantum states of ultracold mat- 
ter [23,[2il,|22. 

In [11] we investigated the ground state of ultracold 
atoms in the optical lattice formed by the interaction 
with the cavity mode. This system combines CQED 



with the many-body physics of strongly correlated ultra- 
cold atoms. In particular, the non-linear dependence of 
the cavity field on the atomic motion opens novel per- 
spectives to the rich scenario of ultracold atomic gases 
in optical lattices. In open space, in fact, these systems 
offer the possibility to reali ze p aradigmatic systems of 
quantum many-body physics 3, [1^ , such as various ver- 



26 . A prominent example is 
27[, exhibiting the superfiuid 



sions of Hubbard models 
the Bose-Hubbard model 
(SF) - Mott insulator (MI) quantum phase transition 
[28| , whose realization with ultracold atoms was proposed 
in Ref. [29] , and demonstrated in Ref. [3^] . In [2J] we ad- 
dressed the question whether and how this transition is 
modified when the atoms are inside a resonator, where 
the optical lattice due to the intracavity-field depends on 
the atomic density. 

In this article we report the details of the derivations 
presented in [2^ and extend them to novel regimes. The 
system we consider consists of ultracold atoms inside a 
resonator, which is driven by a laser. Due to the strong 
coupling between cavity and atomic degrees of freedom, 
the atoms shift the cavity resonance, thus modifying the 
intracavity field intensity. This in turn determines the 
depth of the cavity potential. At ultralow temperatures 
we assume that the atoms occupy only states of the 
lowest band of the periodic optical potential. In this 
regime, we present the detailed derivation of the Bose- 
Hubbard model for atoms in the one-dimensional po- 
tential of an optical resonator, which complements and 
extends the derivation for few atoms by Maschler and 
Ritsch in Ref. [3l|, [13] to large numbers of atoms and 
which is valid in an appropriate thermodynamic limit. 
Using this model we study the SF - MI crossover as a 
function of the system parameters: the chemical poten- 
tial /i, the pump strength and frequency, and the atomic 
density. Assuming the tight-binding regime, we may de- 
scribe the MI states using Wannier functions [s^ , whose 
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form is determined by the optical potential. The Wan- 
nier functions are used then to calculate the coefficients 
in the Bose- Hubbard model, as in standard textbooks. 
We determine the boundaries of the MI states, using the 
strong coupling expansion of Ref . [34] , which is a quite ac- 
curate method for the calculation of the phase-diagram 
of the BH model in ID (25^. 

It must be stressed that the derivation of the Bose- 
Hubbard model in the cavity does involve certain novel 
aspects. Namely, the periodic optical potential depends 
functionally on the atomic density, and hence on the 
Wannier functions. The problem is hence highly non- 
linear: the coefficients of the equations determining the 
Wannier functions depend functionally on the Wannier 
functions themselves, and the latter have thus to be de- 
termined self-consistently. This property has also physi- 
cal implications. In fact, since in our system the coeffi- 
cients of the Hubbard model depend on the atomic den- 
sity, consequently, the diagrams in the plane, where 
t is the tunneling energy, exhibit overlapping, compet- 
ing Mott states, that may even consist of disconnected 
regions for a wide range of parameters. In the vicin- 
ity of the shifted cavity resonance in the strong coupling 
regime, the situation is even more complex: one encoun- 
ters there also dispersive bistable behavior Thus, 
from the quantum optical perspective, this paper inves- 
tigates stability of Mott-like phases, i.e., insulator-like 
states, in an optical resonator. 

This article is organized as follows. In Sec. |TT] we 
present our model that describes a system of two level 
atoms confined along the axis of an optical cavity. In 
Sec. Ill Al we introduce the single-atom Hamiltonian. The 
many-body dynamics including the quantum noise is in- 
troduced in Sec. IIIBI Section III CI is devoted to the 
physical discussion of the role of the various physical 
parameters on the system dynamics. The effective Bose- 
Hubbard Hamiltonian is derived in Sec. lIIII In Sec. lIVI we 
discuss the ground-state properties of our model. We use 
the Gaussian approximation for the Wannier wave func- 
tions, checking carefully its validity. Within the Gaus- 
sian approximation and the strong coupling expansion 
method of H. Monien 34] the stability regions for the 
Mott states are obtained analytically, up to the solution 
of the non-linear self-consistent equation for the width 
of the Wannier functions. Numerical results are reported 
and their physical meaning is discussed in Sec. lIV 5] The 
validity of the approximations is addressed in Sec. IIV CI 
We conclude in Sec.|Vl while in the appendices the details 
of the derivation of the Bose-Hubbard Hamiltonian and 
of the strong coupling expansion method are reported. 



II. THE MODEL 

In this section we generalize the quantum optical model 
of a single atom inside a cavity to the many-body case, 
considering particle collisions at ultralow temperatures 
and quantum statistics. At this purpose, we first intro- 



duce the single-atom dynamics, write then the Hamilto- 
nian in second quantization, introducing atom-atom col- 
lisions, and discuss the basic properties. 



A. Single-particle dynamics 

We consider a single atom of mass M inside a cavity. 
The atomic dipole transition at frequency wq, between 
the ground state \g) and the excited state ]e), couples 
quasi-resonantly with an optical mode of the resonator 
at frequency wave vector k and position-dependent 
coupling strength 

gix) = 5ocos(fca;), 

50 being the vacuum Rabi frequency. The resonator is 
driven by a classical field of amplitude rj and oscillating at 
frequency ujp. We consider the atomic motion along the 
cavity axis, which coincides here with the x-axis, and as- 
sume tight confinement along the radial plane so that the 
transverse motion can be considered frozen out. Atomic 
center-of-mass position and momentum operators are x 
and p, fulfilling the uncertainty relation [x,p\ = iTi. In 
the reference frame oscillating at the frequency LOp of the 
pump field, the normally-ordered Hamiltonian describ- 
ing the coherent dynamics of the atomic and cavity-mode 
states reads 



H 



jc = ^ TiAaCT^ - hl^co^a 



—ihgo cos (x) ((T^a — a^a) — ihi] (a — a^) , (1) 

where = ujp — uja and Ac = LOp ~ lOc are the pump- 
atom and pump-cavity detunings, a (a^) the annihilation 
(creation) operator of a cavity photon at frequency Wc, 
fulfilling the commutation relation [a, a^] = 1, and a = 
](7)(e], fj^ = \e){g\ are the dipole lowering and raising 
operators. Spontaneous emission of the atomic dipole at 
rate 7 and cavity losses at rate k are described within the 
quantum Langevin equations formalism, such that the 
quantum Heisenberg-Langcvin equations for the dipole 
and cavity operators read [l| 



a{t) = (iAc- K)a(t) +.g(x)cr(t) +7] + V2Ka™(t),(2) 
a{t) = [zAa-|]a(t)+g(x)a,a + V7^,r"(t), (3) 

(7^{t) = 2.g(a;) [cj^a -f atfj] - 7 (cr^t) + 1) /2 

+2^{a^r + r^cj), (4) 

where (Jz = cr^cr — cra^ and a*", /™ are the input noise 
operators, whose mean value vanishes and which are 
(5— correlated in time, namely. 



(5) 
(6) 



At large atom-pump detuning the adiabatic elimination 
of the excited atomic state can be performed. Assum- 
ing that the changes in the atomic position are negligible 
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during the time scale in which the atom reaches its inter- 
nal steady state, namely when fc^T ^ we solve 
the Heisenberg-Langevin equations at a fixed valu e of the 
atomic position x. Hence, for |Aa| S> 5o\/ {i^),!, \^c\, 
we set (Jz{t) ~ — 1 in the equations, and obtain tr^ ~ 
ig{x)a^ / Aa- After tracing out the internal degrees of 
freedom, the single-particle Hamiltonian for cavity and 
atomic center-of-mass degrees of freedom reads 



Ho = - — hS, \Uo cos^ (kx) - Ac] a'^a-ihrj (a - a^) , (7) 
2m 



where we have used the explicit form of the cavity spatial 
mode function, and 

= gl/Aa (8) 
is the depth of the single-photon dipole potential. 

B. Many-body dynamics 

We now extend the previous model and derive the cor- 
responding effective Hamiltonian for a gas of N bosons 
at ultralow temperatures. The particle interactions are 
modeled by s-wave scattering. We introduce the field op- 
erators 'i'j{x), with j — g,e labeling the internal 
ground state, such that 



[^,ix),^lix')]^S,, 6ix-x'), 



(9) 



[^,{x),Mx')] = [^]{x),^l{x')] - 0. (10) 

In second quantization the Hamiltonian ([1]) becomes Ti. 
and is decomposed according to H — Ho + Hi, where 

no=J2 fdx^l{x)f^^ + ^u,^]{x)^,ix)\^,ix) 

(11) 

with Uj the strength of the onsite interaction depending 
on the atomic state, and 

Hi = -hAco'a- ihT]{a- a''^) - hAa J dx'^/l{x)'i'eix) 

- ihgo y"dxcos(fca;) [*j;(a;)^'g(a;)a ~ a^^'t(a;)^'e(a;)] . 

(12) 

In the above description we omit to write the Hamil- 
tonian term describing the collisions between atoms in 
different internal states, as we will consider that the ex- 
cited state is essentially empty in the parameter regime 
we choose. The quantum Heisenberg-Langevin equations 



for atomic and field operators read 



-[^g{x),Ho] 



{kx)a^'^e{x) - V7r"^*e(a;) 



«'e(x) 



-50 cos 



--[^e{x),no]+iAa-^e{x) 

h 



(13) 



(14) 



-go cos{kx)^g{x)a - -^^Seix) + y^-^g{x)f 
(iAc - K)a + ri + V2k a™ 



+ga dx cos{kx)-^l{x)'^e{x), 



(15) 



where f™{t) and a"''{t) are the noise operators defined 
in the previous section. Solving the equation for '^e{x, t) 
in the limit of large detuning, | Aa| ^ 7, go\/Jnj, we find 



T / N -90 coa(kx) ^ , , 



(16) 



where the adiabatic approximation lies on the assump- 
tion that ?i| Aa| ^ ksT, as in the single-particle case, and 
we have neglected the input noise term, assuming the de- 
cay rate 7 <C |A|. Substituting this value into Eq. (fT5|l . 
the Heisenberg-Langevin equation for the field is given 

by 



(iAc — K)a + 71 + v2Kaj„ — iJ/oJ^a, 



where 



y 



dx cos^ {kx)'^l ix)^g (x) 



(17) 



is the integral of the density of atoms in the electronic 
ground state, weighted by the cavity spatial-mode func- 
tion squared. We now assume the bad-cavity limit, 
namely the cavity field relaxes to the steady state on a 
much faster time scale than the one in which the den- 
sity of the atomic medium varies. This limit implies 
K ^ ksT/h, and consistency with the previous assump- 
tion imposes |Aq| ^ k ^ kBT/%. In this limit the de- 
pendence of the field on the initial condition is negligible, 
and its solution is essentially the inhomogeneous one that 
can be written as 



a ~ r, F{y). 



(18) 



Here, we have discarded the input-noise terms, as they 
are at higher order in the perturbative expansion and we 
will be dealing with normally-ordered equations, so that 
two-time correlations of the noise operators vanish, see 
Eq. (O. We also introduced the operator 

F{y) = ^- 

^ ' n-i{A,-Uoyy 

(19) 

which is a function of the atom operators in the ground 
state. Substituting Eq. into the equation for the 
ground-state field operator we obtain 



-[^g{x),no]~yC{y,x), 



(20) 



4 



where 

ciy,x) = 7j''Uocos''{kx)FHy)^>g{x)F{y). (21) 

C. Discussion 

Equation ([^D]) shows explicitly the effect of the cou- 
pling with the resonator on the atom dynamics: The 
coupling to the common cavity mode induces a non-linear 
interaction, which enters in the equation through opera- 
tor HH). It is useful to consider the average number of 
photons at steady state riph = (a''^a)st, which we obtain 
from Eq. (fT8)) and reads 

The average number of photons hence depends on the 
atomic density distribution. On the other hand, it deter- 
mines the depth of the confining potential, V « h\Uo\nph, 
and thus the atomic density distribution. In particular, 
the confining potential reaches a maximum for the val- 
ues at which the denominator of Eq. (j22p is minimum. 
From the form of operator one infers that Tiph can 
reach the maximum value when the parameters Ac and 
Uq have the same sign (the operator y is positive valued). 
From Eq. ([8]) this requires that the detunings Ac and A^ 
have equal signs. This property highlights the role of the 
detuning in the dynamics as control parameters. 

We now comment on the parameters required for ac- 
cessing the regime in which the effect of the non-linearity 
will be important, and its consistency with the deriva- 
tion we performed. We first review the important as- 
sumptions, on which our model is based. Spontaneous 
decay is neglected over the typical time scales of the sys- 
tem. This imposes that the effective spontaneous scat- 
tering rate 7', due to off- resonant excitation of the dipole 
transition, fulfills the inequality 7' <^ k. Using that 
7' ~ nphffoT/^a; where riph is the mean value of intra- 
cavity photons, Eq. (|22p . then spontaneous emission can 
be neglected provided that 

2 

"Ph|j|««- (23) 

As our model is based on a single-mode cavity, we also 
require that the detuning between atom and cavity mode 
is smaller than the free spectral range Suj. This reduces 
to the condition 

|A,|«5w (24) 

A further important assumption relies on the relaxation 
time of the cavity field, which has to be much faster than 
the typical time scale of atomic motion. This can be es- 
timated as kbT <C fiK. Finally, in order to have that the 
non-linear effect on riph is sufficiently large for a small 
number of atoms, we have required that Uq ~ k. This 



condition is however not strictly necessary: strong non- 
linear effects can be observed for smaller values of Uq 
when the number of atoms is sufficiently large 

Let us now estimate the number of intracavity pho- 
tons which are usually needed, in order to find the atoms 
in the Mott-insulator state. We consider specifically the 
case in which overlap (bistability) regions between differ- 
ent Mott zone can be observed. In Sec. II VB] we find that 
this occurs at values of the pump amplitude rj ~ 20k. 
This value was evaluated for 50 to 100 atoms in the res- 
onator. Correspondingly, the number of intracavity pho- 
tons is riph ~ 100. From condition we find that 
7 ^ 2|A|a/nph, where we used Uq = k, and which is 
fulfilled for 7 = 27r X 3 MHz and |A|a = 27r x 10 GHz. 
Condition ([24|) is then satisfied when the free spectral 
range Su is of the order of THz. Once |A|a is fixed, 
we find that ~ 2 x 0.1y/K/2TT MHz. Using the value 
K = 27r X 53 MHz [l6|, this requires gg 27r x 700 MHz, 
which is presently at the border of experimental reach. 
However, for smaller values of Uq, say Uq ^ O.Ik, and 
for larger numbers of atoms, say N ~ 1000, the peculiar 
CQED effects on ultracold atoms we predict in this work 
could be well observed for parameter regimes of present 
experiments, see for instance [l6j . 

III. THE BOSE-HUBBARD HAMILTONIAN 

A. Derivation of the Bose-Hubbard Model 

We now derive a Bose-Hubbard type of model for the 
dynamics of the atoms in their self-sustained potential 
when the atoms are well localized in the minima of the 
potential itself. Starting from the assumption that the 
atoms are in a Mott-insulator state, we decompose the 
atomic field operator into the operators b] and bi, which 
create and annihilate, respectively, atoms at the lowest 
band of the potential site centered at x — Xi, according 
to 

^{x)^^w{x-x,)k, (25) 

i 

whereby 'w{x — Xi) are Wannier functions, which are to be 
determined by solving self-consistently the equations of 
motion. The commutation rules of operators bi, b] obey 
the bosonic commutation relations in the regular Bose- 
Hubbard model, where the potential is independent of 
the state of the atoms. We will show that in our case 
this is not apriori warranted, due to the non-linear de- 
pendence of the potential on the atomic density distri- 
bution, which gives rise to a non-linear equation for the 
atomic wave function. However, the bosonic commuta- 
tion relations are still recovered in a properly defined 
thermodynamic limit, which we will identify. 

We rewrite now Eq. ((20)l within this Wannier decom- 
position, 

h = ^[b,,ni''"'^]-iC, (26) 
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where Hq^^^ and C are obtained from Tio and C, re- 
spectively, using the Bose-Hubbard decomposition. They 
read 



BH) 



U 



EoN + EiB + -J2blblbA~fiN, (27) 



and 



C = UorfF^iY) [Joh + Ji{W+i + h-i)] F{Y). (28) 
The couphng matrix elements in Eqs. ([TTll - ipSl) read 



El = I dxw {x - x^)* V^w {x - Xi+i) (29) 



Ji = j dxw {x — Xi)* cos^ [kx] w {x — Xi^i) (30) 
U ^ Ug [dx \w{x)f (31) 



with I = 0, 1 as we keep only on-site and nearest- 
neighbour couplings. In Eq. (j28p we introduced the op- 
erator 



F{Y) 



1 



i(A, - UoY) ■ 



(32) 



where operator Y is the Bose-Hubbard decomposition of 
3^, Eq. (fT7|) . after neglecting couplings beyond the nearest 
neighbors, and takes the form Y = JqN + JiB. 

In order to determine the Bose-Hubbard Hamiltonian, 
we now derive an effective Hamiltonian Hbh such that 
C = [bi,Ti.BH]/^- This is performed in the limit in which 
we can expand operator F in Eq. (|32p in the small quan- 
tity Ji. The details of the derivation are reported in 
App.El The Bose-Hubbard model is recovered for a large 
number of atoms, according to a properly defined ther- 
modynamic limit. We define the thermodynamic limit by 
letting N and the cavity volume to infinity, keeping fi- 
nite the number of atoms per potential site. This implies 
the scaling Uq ~ Additionally, we impose the scal- 

ing ry \/N , which corresponds to keeping the potential 
depth constant as N increases. This scaling corresponds 
to ramping up the pump intensity with ^/N . The Bose- 



Hubbard type of Hamiltonian TicS = Ti-o 
then 



Heff - i^oiV + ^ (n, - 1) - t{N)B + f{N) - iiN , 

(33) 



(BH) 



Hbh is 



where 



t{N) = -El - fnfUoJiF\jQN)F{JoN) (34) 
A,-UoJoN\ 



f{N) = — arctan 



(35) 



We notice that the coefficients of Hamiltonian (|33p are 
operator-valued, hence imposing a Wannier expansion 



such that the coefficients depend on the operator N, 
namely 



w{x - Xi) 



Hence, the commutation relations between the operators 
hi are not the ones of bosonic operators as in the typical 
Bose-Hubbard model. Nevertheless, in the thermody- 
namic limit one finds 



0{l/N). 



(36) 



We therefore perform the Wannier expansion in this 
thermodynamic limit, consistently with the assumptions 
made in order to obtain Hamiltonian (|33p . 



B. The Bose-Hubbard Hamiltonian 

We rescale now Hamiltonian ([55]) in units of the 
strength of the on-site interaction U , which is defined 

in Eq. (pij) . The rescaled Hamiltonian H — H/U reads 



where 



^^ + Eo ^ f{N) 
U NU 



(38) 



contains a rescaled chemical potential, while the tunnel 
parameter 



(39) 



U C/ (k2 + O 

is expressed in terms of the coefficient 

C - A, - UoJoN. (40) 

The higher order terms in JiB, describing long-range 
interaction, have been neglected. Note that the number 

of particles is conserved since [N, H] = 0. We also remark 
that the term f{N)/N tends to a constant and finite 
value in the thermodynamic limit. 

An important physical quantity, which will be useful 
for the following study, is the depth V of the cavity po- 
tential, V — hUoUpYi with riph the number of photons in 
the Bose-Hubbard expansion, Eq. (1^^ . At leading order 
in the expansion in Ji it takes the form 



V 



rfhUo_ 



(41) 



Hamiltonian ([57|) and potential (|¥T|) are the starting 
points of our analysis for the determination of the system 
ground state. 

Let us now make some considerations about the sys- 
tem for a fixed number of atoms A^. From the form of the 
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potential (|4T|) . and in particular from the form of the co- 
efficient C, Eq. we observe that for equal signs of the 
detunings and Ac one can have that ( vanishes. This 
case corresponds to driving the system on resonance, and 
gives a maximum of the cavity mode potential. This res- 
onance situation occurs for atom numbers N that maxi- 
mize the photon number, and gives rise to bistability 3], 
which modifies substantially the properties of the model 
with respect to the standard Bose-Hubbard one. 



IV. DETERMINATION OF THE 
GROUND-STATE 



In this section we determine the ground state of the 
system for a fixed number of particles. Moreover, we 
discuss the situation when the number of atoms is fluc- 
tuating. Our purpose is to identify the parameter regime 
in which the atoms are in the Mott-insulator state. 

Starting from the assumption that the system is in the 
Mott-insulator state, we use the strong coupling expan- 
sion [S^l to verify its validity. In particular, we apply a 
standard degenerate perturbative calculation in the pa- 
rameter t = t/U , and determine the ground-state energy 
Em {no-, t) for the Mott state with uq particles per site, 
and the ground-state energies E± {no, p., t) when one par- 
ticle is added or subtracted to the no-th Mott state. The 
condition 



EM{noJi,i) - E±{no,il,i) = 



(42) 



determines the boundaries fl±{no) of the np-th Mott 
phase as a function of the coupling parameter. For 
/x+(no) > fi-{no) the region between the two chemical 
potentials determines the Mott zone. The Mott state 
gets unstable as the parameters are varied such that 
/2+(no) = fl-{no) and finally /i+(no) < M-(^*o)- In this 
section we determine the boundaries of the Mott state in 
a diagram, in which we plot /i as a function of relevant pa- 
rameters. We remark that, in the typical Bose-Hubbard 
model, when the system exits the Mott phase, then it is 
in a superfluid state. In our case, this is probably verified 
in most cases, which we will discuss individually. 

Finally, the parameter t can be controlled by vary- 
ing the pump amplitude ?7, which is straightforwardly 
related to the number of of photons inside the cavity 
and hence to the height of the potential. Alternatively 
it can be changed by varying the atom-pump detuning 
Aa and the cavity-pump detuning Ac, which enter in the 
dynamics through the coefficient (j40p in the denomina- 
tor of Eq. and correspond to changing the refractive 
index of the atomic medium. 

In the following wc first study the functional depen- 
dence of the integrals on the system parameters using 
the Gaussian ansatz. We then determine numerically the 
regions of the Mott-insulator state in the diagram where 
the chemical potential /2 is studied as a function of the 
pump intensity rj. 



A. Coefficients in the Gaussian approximation 

We determine the boundaries of the Mott-insulator re- 
gions using the Gaussian approximation, hence replacing 
the Wannier functions by Gaussian functions in the inte- 
grals (p9| - (pT|) . Thus, the Wannier functions are replaced 
with Gaussian functions such that 



w{x - Xi) K. wg{x - Xi) = (tto-^) 



^1/4, 



(43) 



where a is the width to be determined. This treatment 
allows us to identify the dependence of the coefficients on 
the physical parameters, reproducing with good approx- 
imation the results obtained with the Wannier functions 
in the parameter regimes we discuss in Sec. IIV CI In 
particular, we modify the Gaussian functions in order to 
fulfill the orthogonality condition, 



da; 'w'q{x — Xi)wQ{x — xj) — Si 



In this way we avoid small, but unphysical contribu- 
tions. Let K be the number of lattice sites. The width 
a of the Gaussian functions is found from the depth V 
of the cavity-mode potential, Eq. (|4ip . In particular, 
— hj ^2m\y\k. In order to determine the boundaries 
of the Mott states in the diagram of /i as a function of 
r], we determine the coefficients for the three cases (1) 
iV = Kno + 1, (2) TV = Kuq and (3) N = Kno - 1, 
and introduce the subscript (i) with i = 1,2,3 for the 
corresponding coefficient. We evaluate the integrals in 
Eqs. (|29|) - ([3T|l for these three cases and express them as 
a function of the dimensionless parameter 



yi^) = = ^ErI\V^,)\, (44) 

where Eji is the recoil energy. In term of , they read 
Er 



Eo(i} — 



Ei(i) = 
j± - 



[l - sign(Aa)exp (-?/(.,))], 



(45) 



(46) 



4 

sign(Aa 



exp 



■ exp 



4:ERa, 



/27rA„ 



4) (47, 

yw), (48) 

(49) 



where is the scattering length, Aj^^ is the atomic wave 
packet transverse width and the sign ± depends on the 
sign of A(j. In the limit J^.^ 3> l^j^j-,! the potential 
amplitude according to (j4ip is given by 



V, 



(i) 



ri^nUo 



(50) 
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As J^.-j depends on V(i) which, on the other hand, de- 
pends itself on J^^ , the above equations must be solved 
self-consistently. This is a consequence of the atom- 
density dependence on the coupling parameters. In par- 
ticular, for Aq > (atoms at the nodes), Jo(i) — *■ in 
the strong pumping limit, rj ~* oo, and the results be- 
come independent of the number of atoms. On the other 
hand, if < (atoms at the antinodes) the parameter 
Joii) ~^ 1 for sufficiently large pumping, and the non- 
linearity is strongest. 

Within this treatment we determine the nearest- 
neighbor coupling parameter, which is given by 

= ^^r^'"^ 2^"^'") ' (51) 

where U = 2h^ agk / (\/2TTmAyz) . For rj ^ oo the poten- 
tial |V(i)| oo and consequently y^j) — s- 0, and hence 

The perturbative calculation of the boundaries fl±{nQ) 
is sketched in the App. |B] At third order the result reads 

^(12) 

jj.+ {no) = no H —Kno{no - 1) - i(i)2(no + 1) 

+ tl)nl - - *(2)^) -^Knoino + 1) 
+ ifi)no(no + l)(no + 2), (52) 
fl-{no) = (no - 1) —Kno{no - 1) + i(3)2no 

- i(3)K + 1)' + (^<^3) - ^(2)^) ^Knoim + 1) 
+ if3)no(no + l)(no-l). (53) 
Here, [/(^a) = 1 - U(2)/U(i), fi+{no) = /i+(no)/C/(i) and 

B. Numerical results 

In this section we study the regions of the Mott- 
insulator state as a function of the chemical potential 
and of the inverse pump amplitude 77^^. The boundaries 
are determined by numerical evaluation of Eqs. (|29p - (PT|) 
using the modified Gaussian functions. The atomic pa- 
rameters we choose correspond to *^Rb atoms with scat- 
tering length Os = 5.77 nm and atomic transition wave- 
length A = 830 nm. The optical potential has K lattice 
size and the transverse width of the atomic wave packet 
is Aj^ — Az = \/ Ayz — 30 nm. We evaluate the "phase 
diagrams" for K = 50 — 10000 at fixed number of atoms 
N, scaling TV so to keep the atomic density constant. The 
results for the Mott zones agree over the whole range of 
values, so in the figures we report the ones obtained for 
-ftT = 50 for different values of the detunings. 

Figure [T] displays the first four Mott zones for (a) 
Aa < and Ac = k and (b) Aq > and Ac = 0, as 




0.01 0.02 0.03 0.04 0.05 0.06 0.07 0.08 
K/T| 



(b) 41 

2 

n 















1 






A =0 

u;=5k 





0.05 0.1 0.15 



kJt\ 

FIG. 1: Phase diagram, showing the Mott zones in the scaled 
chemical potential-cavity pump plane jl-ri~^ . The parameters 
are (a) Uo — —5k (Aa < 0, atoms at the antinodes), Ac = 2k 
and (b) Uo = Sac (Aa > 0, atoms at the nodes), Ac = 0. The 
overlap and reappearance of the Mott zones originate from 
the non-linearity of the system. The dotted lines correspond 
to the boundaries of the covered zones. 



a function of the dimensionless parameter k/t]. Interest- 
ingly, the extension of the Mott zones seems to decrease 
roughly as Hq^ in both cases. We first analyse the case 
displayed in Fig.[lja). For A^ < the atoms are trapped 
at the maxima of the intracavity field. Hence, the cou- 
pling with the cavity mode is maximum when the con- 
finement is very tight. Here, for large values of the pump 
intensity (i.e., for small values of k/tj) the Mott zones at 
different values of no show some overlap. This overlap 
is a cavity QED effect, in fact Mott states with larger 
number of atoms per site are favored as they increase the 
coupling strength to the cavity mode, and thus the depth 
of the potential. The overlap is only at the border of the 
boundaries, as atom-atom collisions compete with this 
effect. In Fig. [Ijb) the detuning A^ > 0, and the atoms 
are hence trapped at the nodes (the zeroes) of the intra- 
cavity field. Hence, the coupling with the cavity mode is 
minimum when 77 ^ oo. Indeed, here we observe that for 
large values of 77 (small values of k/t]), the Mott zones 
almost do not overlap. However, for smaller values of 
77 they exhibit an "exotic" behavior: overlap, disappear 
and reappear. 

Further insight is gained in Fig. [2l where we study 
the depth of the cavity potential as a function of the 
pump parameters. The curves displayed in Fig. [2Ka) 
correspond to the parameters of the phase diagram in 
Fig. [ija). Here, one observes that the potential ampli- 
tude increases monotonically as y ~ 77^ in the parameter 
regime where the non-linearity is weak. Correspondingly, 
the width crjj) of the Wannier functions, giving atomic 
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localization at the minima of the potential, decreases 
smoothly as a^i) cx iVf^)!"^/"* ^ l/y^, see Eq. (j44)) . For 
larger pump strengths, when the non-linearity becomes 
important, the behaviour is slightly changed. The curves 
in Fig. IHb) correspond to the parameters of the phase 
diagram in Fig. [Ijb). Here, one finds that the potential 
depth increases rapidly where the corresponding Mott 
zones exhibit a minimum in the value of fl. Correspond- 
ingly, the width cr(j) c>c |V(i) |~^^* diminishes rapidly. This 
behaviour changes at the value of 77 where the potential 
gradient increases abruptly. In this regime ai varies very 
slowly. This can be understood as a competition be- 
tween the cavity field, which tends to localize the atoms 
at the minima, and the atomic quantum fluctuations: 
When the potential is sufficiently high to trap the atoms 
within a small fraction of the wavelength, the cavity field 
is pumped more effectively. 




0.05 0,1 0.15 0.2 



FIG. 2: (Colour online) The potential amplitude |V(i)|, in 
units of Er and in log-scale, as a function of the inverse pump- 
ing K/rj, where the curves in (a) and (b) have been evaluated 
in the parameter regimes of FigslTJa) and[TJb), respectively. 
When the non-linearity is strong the potential amplitudes dif- 
fer from the linear situation where \V\ ^ rf . The average 
number of cavity photons is found by multiplying the rescaled 
potential depths in the plots by the factor /„ = Er/\hUo\, 
which here is /„ « 0.006. 



We now consider the situation in which the detunings 
Aa and Ac have the same sign. In this case the pa- 
rameter C{N) in Eq. vanishes when the condition 
Jg*^.) = Aa/UoN is fulfilled, whereby < J+^j < 1/2 and 

1/2 < Jp"^.^ < I, see Eq. This resonance condition 

gives rise to bistability, leading to an abrupt change of 
the potential depth. As a consequence, the Mott state 
may become unstable. The upper plot of Fig. [3] displays 



the potential amplitudes V(i) /Er for one atom per site as 
function of tj/k for Uq = —k and Ac = —45k, exhibiting 
the typical functional behavior of bistability. The lower 
plot displays the corresponding phase diagram. For the 
case of one atom per site, the Mott ground-state will 
suddenly disappear for rj ^ n when the pumping is adia- 
batically lowered. Clearly the first "jump" occurs in the 
potential V(i) (corresponding to the lowest atom density), 
and comparing the two plots one finds that this takes 
place exactly when the first Mott zone suddenly ends. 
The system most likely jumps into a state where higher 
Bloch bands are populated. In this case the single-band 
and Gaussian approximations break down. 




0.2 0.4 0.6 0.8 1 1.2 
k/t] 



FIG. 3: The bistability behaviour of the potential amplitudes 
V(i) as a function of n/ri (upper plot) and corresponding phase 
diagram fl — 77"^ (lower plot) for Ac = —45k, Uo ~ — k. At 
no = 1 the Mott region suddenly ends for ^ k, where the 
corresponding potential V(i) jumps to a lower value. Here, the 
system most likely is in a state where higher Bloch bands are 
populated, due to non-adiabatic effects and the small poten- 
tial depth of this solution. 

The overlapping of the Mott zones and the bistability, 
which we observe in the phase diagram, are novel features 
when compared with the typical scenario of cold atomic 
gases trapped by an external potential. Let us first dis- 
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cuss on the existence and uniqueness of the ground-state. 
When the Mott-insulator state is stable, given the num- 
ber of atoms N, the ground state is fully determined once 
the atomic density p = N/K is fixed. Outside of the Mott 
zones we expect superfluidity in the parameter regimes in 
which there cannot be optical bistability (detunings with 
opposite signs). In the situation of multiple solutions of 
the Eqs. the system will most likely be found 

in the one solution which minimizes the energy. 




FIG. 4: Extended phase diagram of fig. [T] (a), where the 
atomic density p has been included as a third axis. Here, 
the contour lines correspond to a fixed atom density p, such 
that for a given p the scaled chemical potential depends on 
the pumping strength r) according to this particular contour 
line. The dotted contour curves indicate the lines with exactly 
no atoms per site. The projection of the Mott-zones onto the 
p. — rj~^ plane is shown. 



A more complete picture of the phase diagram can be 
obtained by considering the dependence on the atomic 
density. The strong coupling method for higher orders is 
cumbersome once the number of added/subtracted par- 
ticles to the Mott states becomes larger than one. How- 
ever, the first-order corrections are still easily obtainable 
for any atom number. In fig. 3] we present schemati- 
cally the extended phase diagram of fig. [1] (a) where the 
atomic density has been included as a third axis. This 
diagram has been obtained by fitting the intermediate 
lines between the Mott zones, and by verifying that it 
reproduces the first order calculations for small values of 
k/t]. We remark that the "overhangs", corresponding to 
the overlapping zones in Fig. [^a), constitute the novel 
feature, which we encounter in this model as compared 
to the standard Bose-Hubbard model. 

When the number of atoms is not fixed [35] , the atomic 
density may take multiple values where the phase dia- 
gram exhibits "overhangs". For sufficiently long times, 
we expect that the system will be found in the number 
of atoms such that the energy is minimized. This implies 
also that there may be a competition between a Mott 
and a superfluid state at two different values of the den- 
sity, which happen to be at similar energies. Keeping 



this situation in mind, we restrict our analysis to differ- 
ent and overlapping Mott states, and compare their en- 
ergy. Figure [5] displays a phase diagram on the p — ri~^ 
plane, whereby the Mott states with higher energy are 
plotted on top of the ones with lower energy. We observe 
that for large pumping strengths the Mott states with 
a higher number of atoms no have in general a greater 
energy, while for lower or moderate pumping strengths 
this is not necessary true. For example, the end of the 
Mott zone with uq — 4 has smaller energy than the cor- 
responding one for the uq — 3. This is a pure CQED 
effect. 




FIG. 5: Phase diagram on the p — rj~^ plane, reporting the 
first four Mott zones. The Mott states with higher energy are 
plotted on top of the ones with lower energy. Typically, for 
large pumping the Mott zones with a large number of atoms 
no per site have the highest energy. For moderate pumping 
this is not necessarily true as seen by comparing for example 
the third and fourth Mott zones at around n/r] ~ 0.05. The 
relevant parameters are reported in the inset. 



C. Validity of the approximations 

We now discuss the regime of validity of the calcula- 
tions, from which we extracted the phase diagrams pre- 
sented in this section. The derivation of the system cou- 
pling parameters relies on the assumption Jq ^ I'^i^l- 
The maximum value of the ratio \ J^\/ Jq ~ 0.056 occurs 
at \V\ « Er/2, hence the nearest-neighbor coupling is at 
least 17 times smaller than the on site coupling. The ex- 
pansion to first order in JiB of Eq. (j32p is motivated for 
any number of atoms since the perturbative parameter 
A = JiB/JqN ^ Ji/Jo is strictly smaller than unity. 

The values of the chemical potential, as in Eq. (|52|). 
are derived from a third-order perturbation expansion 
in the parameter t = t/U and it is expected to break 
down for large t. We verified that in general t < 0.25. 
Moreover, we compared the phase diagrams with the ones 
obtained by truncating at the second order in t, and could 
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verify that they do not differ substantially one from the 
others. We remark that the perturbation calculations 
are carried out assuming periodic boundary conditions, 
while the system here studied has a fixed number of sites, 
K — 50. We checked the validity of the assumption by 
comparing the results obtained for different lattice sites, 
up to 10000, keeping the density fixed. 

As it concerns the tight-binding approximation (i.e., 
only including nearest neighbor couplings), the single- 
band approximation (i.e., expanding the field operators 
^'(f) and ^^x') using only the lowest band Wannier 
functions), these are both related to the regime of valid- 
ity of the Gaussian approximation. Within this approx- 



3> 1, also 



imation one finds |Ji/J„| = exp 

indicating validity of the tight-binding approximation in 
this regime. Figure[BJa) displays the difference A—Atba 
between the width A of the first Bloch band, obtained 
from diagonalization of the single-particle Hamiltonian in 
Eq. ([7]), and the width Atba evaluated with the tight- 
binding approximation, as a function of y'^- FigurelHI^b) 
displays the difference between the coupling parameters 
obtained by using the Wannier functions and the mod- 
ified Gaussian functions as a function of y~^- We note 
that for values y^^ < 1 the validity of both the tight bind- 
ing and Gaussian approximation visibly breaks down. 
This has also been verified by recalculating some of the 
above phase-diagram using the Wannier functions. 




2 -1 3 

y 



FIG. 6: Check of the validity of the involved approximations. 
The upper figure shows the difference A — Atba, in units 
of Er, as a function of , where A is the width of the 
first energy band and Atba ~ ^tw, with tw the coupling 
element obtained from the corresponding Wannier functions 
in the tight-binding approximation. The lower figure displays 
tw —to as a function of y~^, where to is the coupling element 
given by the Gaussian approximation. The quantities are 
numerically derived from Hamiltonian (0, where, in scaled 
units, there is only a single parameter of interest, namely the 
dimensionless potential amplitude or equivalently y. 



V. CONCLUSIONS 

We have shown that ultracold bosonic atoms inside a 
resonator may form stable insulator-like states, and thus 
enter the Mott phase, which is sustained by and sustains 
the cavity potential. The low temperature properties of 
the system are determined by the competition between 
the quantum electrodynamic effects and the quantum 
fluctuations of the atomic matter waves. This compe- 
tition gives rise to a non-trivial dependence of the re- 
gions of stability and of the collective atomic states on 
the system parameters. Since the cavity potential de- 
pends on the state of the atoms, the behavior of the 
ultracold atomic gas in the cavity differs hence signifi- 
cantly from the one which is encountered in open space. 
We have derived the Bose-Hubbard Hamiltonian for the 
cavity confined system, and have shown that the coeffi- 
cients of this Hamiltonian depend explicitly on the num- 
ber of atoms. We have determined regions of parameters 
where the atomic insulator states are stable, predicting 
the existence of overlapping stability regions for compet- 
ing Mott states. Bistable behavior is encountered in the 
vicinity of the shifted cavity resonance, controlled by the 
pump parameters. 

Our theory allows us to determine the state of the 
atoms when their number is fixed, while for fluctuat- 
ing, non-fixed atom number, in general, the system will 
choose the state of minimum energy. This will also take 
place when an external inhomogeneous potential, such 
as a harmonic trap potential, is additionally applied to 
the atoms. In such case we envision a possibility of hys- 
teresis effects in the harmonic potential, when the fre- 
quency of this potential is slowly increased and, subse- 
quently, slowly decreased. However, the question how the 
presence of an inhomogeneous potential will modify the 
insulator-like states requires further careful studies, since 
the state of the system depends in a highly nontrivial way 
on global parameters, which in turn determine the local 
density of the atoms and the intracavity field. The con- 
dition that atoms may affect locally the potential, hence 
giving rise to phonon-like feature [? ] , may be reached in 
multi-mode resonators, allowing for localized polaritonic 
excitations [13, . Further novel features are expected 
when fermions are considered instead of bosons. These 
questions will be tackled in future works. 
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APPENDIX A: DERIVATION OF THE 
EFFECTIVE HAMILTONIAN 



We consider Eq. 



and rewrite it as 



bi = T-[be,no]-iC, 
in 



(Al) 



where C is defined in Eq. ^ and Y = JqN + JiB. 
We aim at finding an effective Hamiltonian TisH of the 
Bose-Hubbard form, such that 

c^[be,nBH]/h 

in some thermodynamic Umit to be identified. 

We expand now operator C at first order in Ji , assum- 
ing Ji <C Jo, as it is verified in the Mott-insulator state, 
using [N,B] — and 

F{JoN + JiB) « F{JoN) + JiBF'iJoN), (A2) 
where we have introduced the notation 



At first order in Ji, we find 



(A3) 



C = Uor^^F^{JoN)[Joh + Mbi+,+b,_,)]FiJoN) 
+UoV^F'\JoN)JiBJobiF{JoN) 
+UoTfF\JoN)Jnbi.hBF'{JoN) + 0{JI). 

Let us now consider the commutation relations between 
the various operators entering this expression. We note 
that 

[b,,F{JaN)] = {F{Jo{N + l))-F{JoN))bi 

= F'{JoN)Jobi + 0{l/N^) (A4) 

and it is hence of order 1/N. Similarly, the commutator 
[bi, B] = 6^+1 + bi-i is at higher order in the expansion 
in 1/N. Henceforth, we can rewrite 

C = Uo7j^F\JoN)[Jobe + Jiibi+i+bi^i)]F{JoN) 
+Uov^F'^{JoN)JobeJ,BF{JoN) 
+UoV^F^iJoN)JiBJobeF\JoN) 
= \bi,T-LBH], 

where 

Ubh = nifUoJiF\jQN)BF{JoN) + G{JaN) (A5) 

and operator G{JqN) has to be determined from the 
equation 

[bi, G{JqN)] + UorfF^{JoN)JobiF{JoN) = 0, (A6) 



which is valid at the considered order in the expansion in 
1/A^. At leading order in 1/iV, Eq. (KE\i is a differential 
equation, such that G'iJoN) = -Uqv'^F'i {JoN)F{JoN). 
Using the explicit form of operator F{x), Eq. (|32p . we 
find 

G'(x) = -UovVi^^ + (Ac - Uox)^) 
which gives 

G(.x) = ^arctan('^^^') (A7) 



and finally the effective Hamiltonian in Eq. ([33 



APPENDIX B: PERTURB ATIVE DERIVATION 
OF THE ZONE BOUNDARIES 

We consider Hamiltonian 



H = -i{N)i3 



U{N) 



-l)-/iiV, (Bl) 



as given in Eq. (|37p . This Hamiltonian differs from the 
standard Bose-Hubbard Hamiltonian, as the coefficients 
depend on the operator N. We apply now to Eq. (jBl[) 
the method of Ref . [3j| , which allows to determine the re- 
gion of stability of the Mott-insulator states. The method 
consists in a perturbative expansion in the parameter 
t, which is assumed to be small within the parameter 
regime of interest. In this limit, for large onsite interac- 
tion strength U (hard core limit), in the optical lattice 
the configuration which is energetically favorable has the 
smallest number of atoms per site. For a lattice of K 
sites and N = Kuq + j atoms, with j < K, there will 
be either no or tiq 1 atoms per site. Clearly, when 
N = Kuq atoms {j = 0), there exists only one possible 
ground-state, while for N > Kuq several ground state 
configurations exists, and one has to apply degenerate 
perturbation theory. 

The ground-state of Hamiltonian (jBl[) is found after 
imposing periodic boundary conditions, and diagonal- 
izing operator B in the momentum representation. At 
t = the ground-state for N = Kuq is given by 



|*o("-o)) = \nQ,no, ...,no), 



(B2) 



corresponding to tiq atoms per site, while for N — Khq- 
j, with j > 0, they are defined by the relation 



l*.K)>=4,|vI/,_i(no)) 



where 



1 ^ 



Ank-ia n 



\/nn + 1 



(B3) 



(B4) 



creates one particle in a site starting from the low- 
est energy states. There is an analogous state for one 
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hole. Here, a = Tr/fc is the distance between neigh- 
boring sites, and the wave vector kj = 2'Kj/Ka, with 
j = — -y,— -y + l,...,-y — 1 (assuming K even for simphc- 
ity). The ground state energy is calculated applying per- 
turbation theory in third-order in t to this unperturbed 
basis. Due to symmetry, only zeroth and second-order in 
the perturbation of t{N)B contribute to the ground-state 
energies of the Mott-insulator state. For TV = Kno one 
finds 



^(2) 



(B5) 



(2) 



2Knoino + 1) 



while the SF energies for the added particle/hole energies 



E+{no) = %i[ifno(no-l) + 2no]-Ai(i)(i^"o + l) 



Here we have used the subscript {i) corresponding to the 
three different cases, N = Kuq and N = Kuq ± 1, see 
Sec. lIVI The limit of stability of the Mott-insulator state 
is found when the states l^'M(fio)) and |4'±(no)) are de- 
generate. The conditions Em (no) — E^{nQ) — and 
i?j\/(no) — E^Itiq) ~ determine the boundaries /i±(no) 
of the Mott states in the phase diagram fl — t, thus ob- 
taining the results in Eqs. 



-t(i)2(no + l) 



"(1) 



[2Kno{no + 1) 



-jj^na{no + l)(na + 2), 
^(1) 



-t(3)2no - ^ [2Kno{no + 1) - (no + 1)^] 



••(3) 



no(no -I- l)(no - 1). 



(B6) 
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